% This function generates power curves for LOCA, CCE, IM, and CRS.
function [output] = output_power_curve(setting,n_sim,alpha_sig)

cd ../result
p_adj_vec = zeros(n_sim,1);
l_power_curve = 50;
p_power_curve_mat = zeros(n_sim,l_power_curve);
p_null_vec = zeros(n_sim,1);

i_real = 0;

for i_sim = 1 : n_sim
    s = ['results',num2str(i_sim),'.mat'];
    if exist(s, 'file') == 2
        i_real = i_real+1;
        load(s)
        eval(['output = output_',setting,';']);
        p_power_curve_mat(i_real,:) = output.p_power_curve;
        p_adj_vec(i_real) = output.p_adj;
        p_null_vec(i_real) = output.p_null;
    end
end
cd ../output
p_power_curve_mat(i_real+1:end,:) = [];
p_adj_vec(i_real+1:end,:) = [];
p_null_vec(i_real+1:end,:) = [];

power_curve = mean(p_power_curve_mat<repmat(p_adj_vec,1,l_power_curve));

% adjused p-value threshold such that empirical rejection rate is alpha_sig
alpha_adjusted = quantile(p_null_vec,alpha_sig); 
% adjusted power
power_curve_adjusted = mean(p_power_curve_mat<alpha_adjusted+eps);


%% OUTPUT
output.power_curve = power_curve;
output.power_curve_adjusted = power_curve_adjusted;


